#encoding=utf-8

# 已知椭圆曲线加密Ep(a,b)参数为
# p = 15424654874903，a = 16546484，b = 4548674875，
# G(6478678675,5636379357093)，私钥为k = 546768，求公钥K(x,y)

# 扩展欧几里得算法是欧几里得算法（又叫辗转相除法）的扩展。通常谈到最大公因子时, 我们都会提到一个非常基本的事实: 
# 给予二整数 a 与 b, 必存在有整数 x 与 y 使得ax + by = gcd(a,b)。因此，有两个数a,b，对它们进行辗转相除法，可得它们的最大公约数；
# 然后，收集辗转相除法中产生的式子，倒回去，可以得到ax+by=gcd(a,b)的一组整数特解。

def exgcd(a, b):
    if b == 0:
        return 1, 0, a
    else:
        x, y, q = exgcd(b, a % b)#递归直至余数等于0(需多递归一层用来判断)    
        x, y = y, (x - (a // b) * y) #辗转相除法反向推导每层a、b的因子使得gcd(a,b)=ax+by成立     
        return x, y, q

# 扩展欧几里得求逆元

def invert(a, p):
    x, y, q = exgcd(a, p)
    if q != 1:
        raise Exception("No solution.")
    else:
        return (x + p) % p  # 防止负数

#扩展Eulid算法（求乘法逆元）
def ExtendedEulid(a:int,b:int):#ax=1modb,得到b在模a下的乘法逆元
    def ExtendedEulid0(a:int,b:int):
        if b==0:#边界条件
            return 1,0,a
        else:
            x,y,gcd=ExtendedEulid0(b,a%b)#递归
            x,y=y,(x-(a//b)*y)#递推关系，左端为上层
            return x,y,gcd#返回第一层的计算结果。
        #最终返回的y值即为b在模a下的乘法逆元
        #若y为负数，则y+a为相应的正数逆元
    n=ExtendedEulid0(a,b)
    if n[1]<0:
        return n[1]+a
    else:
        return n[1]



def ECCnegative(x, y, p):
    return x, p-y


def ECCadd(x1, y1, x2, y2, p, a):
    if(x1 == x2):
        if((y1+y2) % p == 0):
            return -1, -1  # infinity return
        elif(y1 != y2):
            print("Input Error, Please check the input!")
            return -2, -2  # error input return
        else:
            k = ((3 * x1 * x1 + a) * invert(2 * y1, p)) % p
            x3 = (k * k - x1 - x2) % p
            y3 = (k * (x1 - x3) - y1) % p
            return int(x3), int(y3)
    else:
        k = ((y2 - y1) * invert(x2 - x1, p)) % p
        x3 = (k * k - x1 - x2) % p
        y3 = (k * (x1 - x3) - y1) % p
        return int(x3), int(y3)


def ECCmul(k, x, y, p, a):
    tempx, tempy = -1, -1
    nowx, nowy = x, y
    while(k != 0):
        mark = k & 1
        k = k >> 1
        if(mark == 1):
            if(tempx == -1):
                tempx, tempy = nowx, nowy
            else:
                tempx, tempy = ECCadd(tempx, tempy, nowx, nowy, p, a)
        nowx, nowy = ECCadd(nowx, nowy, nowx, nowy, p, a)
    return tempx, tempy


k = 546768
p = 15424654874903
a = 16546484
b = 4548674875
print(ECCmul(k, 6478678675, 5636379357093, p, a))
